Dissipative Dynamics with Trapping in Dimers 
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The trapping of excitations in systems coupled to an environment allows to study the quantum to classi- 
cal crossover by different means. We show how to combine the phenomenological description by a non- 
hermitian Liouville-von Neumann Equation (LvNE) approach with the numerically exact path integral Monte- 
Carlo (PIMC) method, and exemplify our results for a system of two coupled two-level systems. By varying 
the strength of the coupling to the environment we are able to estimate the parameter range in which the LvNE 
approach yields satisfactory results. Moreover, by matching the PIMC results with the LvNE calculations we 
have a powerful tool to extrapolate the numerically exact PIMC method to long times. 

PACS numbers: 05.60.Gg, 05.60.Cd, 7I.35.-y 



I. INTRODUCTION 

Recent years have seen growing interest in coherent energy 
transfer. For instance, it has been pointed out that photosyn- 
thesis might benefit from quantum mechanical features of the 
transfer of excitations created by the incoming solar energy 
d]. A series of papers has modelled coherent dynamics in 
the light-harvesting system of the photosynthetic unit show- 
ing that the coupling to an environment does not necessarily 
destroy all coherent features - even at room temperature - but 
also can support the coherent transfer of excitations [2J] . The 
majority of these studies use the Lindblad form of the Liou- 
ville von Neumann equation (LvNE) for the reduced density 
operator of the system where the environmental degrees of 
freedom have been traced out [3]. However, this approach is 
only valid in a limited parameter range of the coupling to the 
environment. 

In contrast, rapid experimental advances allow to manipu- 
late and control ultra-cold atoms to a large extent. This offers 
the possibility to study coherent transport and the effect of en- 
vironment changes (e.g., an increase in the temperature). An 
ideal system to study the dynamics of excitations are (frozen) 
Rydberg gases |4J], whose atoms can have complex spatial 
arrangements, for which the coherent dynamics can be effi- 
ciently modelled by continuous-time quantum walks J5D. 

Moreover, it is possible to adjust the properties of specific 
("special") single atoms such that the excitation gets to be ab- 
sorbed at these atoms 10,01. In this way they mimic the reac- 
tion center (RC) of photosynthesis, where the excitation gets 
trapped and further processed. Both systems can be viewed 
as being donor-acceptor units, where the excitation is created 
at the donor ("normal" Rydberg atoms, light-harvesting com- 
plex) and gets absorbed at the acceptor ("special" Rydberg 
atoms, RC). The decay of the excitation at the acceptor allows 
to globally monitor the transport dynamics. 

While the LvNE approach allows for a phenomenological 
modelling, other methods treat the system and the coupling 
to the environment in a microscopic manner. Our method 
of choice is the path integral Monte Carlo (PIMC) technique, 
see, e.g., fil [T^l . which can be applied for arbitrary system- 
environment coupling strengths. However, unlike the LvNE 
approach, the real-time PIMC method is plagued by the noto- 
rious dynamical sign problem [14], which significantly limits 



the time scales accessible to PIMC simulations. However, by 
combining the LvNE and the PIMC methods we are able to 
study excitation dynamics and trapping over large timescales 
and improve the numerical accuracy of our results. 

Our model system is a dimer, represented by two coupled 
two-level systems (TLS), one of which acting as trap, the 
dimer being coupled to the environment. By assuming not 
too strong couplings to the environment and a single initial 
excitation of one of the TLS, one can map the two TLS onto a 
single TLS, if, without the trap, the probability of finding the 
excitation in the system is conserved [9]. We note that various 
systems with, e.g., radial symmetry and a trap in the center can 
effectively be mapped onto the dimer if initially the excitation 
is homogeneously distributed over the periphery ifloll . 



II. COHERENT EXCITON TRAPPING 

In general, we consider the Hamiltonian of a network of 
N nodes, where each node represents a single two-level sys- 
tem. Let Hq be the Hamiltonian without traps. The acces- 
sible Hilbert space is completely characterized by the basis 
states \j) which are associated with excitations localized at 
the nodes j ~ 1, . . . , N. Within a phenomenological ap- 
proach, the Hamiltonian, which incorporates trapping of ex- 
citations at the nodes m e Ai, M. C {1, . . . , N}, is given 
by H = Hq — iT, where iF = iTj^m l m )( m l is the trap- 
ping operator, see Ref. [6] for details. As a result, H is non- 
hermitian and has N complex eigenvalues, Ei = e; — i<ji 
{I = 1 , . . . , N) where 7/ > 0, and N right and N left eigen- 
states, denoted by |<£>/) and respectively. The transition 
probability from node j to node k is then given by 



E 



exp(- 7/ i) exp(-ieit)(k\$i)($i \j) 



(1) 



so that the negative imaginary parts 7/ of Ei determine the 
temporal decay. The mean survival probability Ii(t) of an ex- 
citation in the presence of M trap nodes is a global property 
of the network and is defined as 



n(t) = 



1 



N — M 



(2) 
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i.e., 11(f) is the average of iTkj{t) over all initial nodes j and 
all final nodes k, neither of them being a trap node. 



eigenvalues, i.e., 7+ = 7_ = T/2. The bi-orthonormalized 
eigenstates of H are of the form 



A. Liouville-von Neumann equation 



1 



^J2 cos < 



e ±i<t>/2 



(7) 



The Schrodinger equation can be recast into the Liouville- 
von Neumann equation (LvNE) when considering the density 
operator p. For Hermitian Hamiltonians Hq one has p = 
—i [Ho, p\ , where [•, •] is the commutator. Now, for the non- 
hermitian Hamiltonian H one obtains 



and 



P = -i[H ,p] - {T,p}, 



(3) 



where {•, •} is the anti-commutator. 

Introducing the coupling to the environment naturally com- 
plicates the situation. However, under certain conditions one 
can employ the so-called Lindblad form of the LvNE, where 
the Lindblad operators specify the coupling fl. Consider now 
Lindblad operators which can be written as y/XjLj, where Xj 
is a fixed decay rate. Assuming all rates to be equal, Xj = A 
for all j, the LvNE reads 

P=- 1 [H ,p}-{T, P }-2XJ2(p-(j\p\j))l j . (4) 



The rate A can be estimated from the spectral density 
J(w) describing the environment within the Caldeira-Leggett 
model [12] at a given temperature T. Taking J(u>) = 
2-kcluj exp(— ui/lu c ) and using the Markov approximation one 
arrives at A = iraksT yfl. One has to bear in mind that 
Eq. is an approximation with a limited range of validity: 
For a very large coupling strength A, Eq. (0]l leads to the quan- 
tum Zeno limit rather than to a classical master/rate equation. 
In the following we will consider Lindblad operators which 
are given by projection operators of the type Lj = \j) (j\ Jit]. 



III. DIMERS WITH TRAPS 

In the sequel we will consider a dimer which is coupled to 
an external bath. This system allows to solve Eq. analyt- 
ically and, moreover, to compare the approximate LvNE re- 
sults to the numerically exact PIMC calculations. The Hamil- 
tonian of the dimer without any coupling to the surroundings 
can be expressed through the Pauli matrices a z and <x x , 



H = E 1 - Vcr r 



(5) 



where E is the onsite energy, which we choose to be the same 
for both nodes, and V is the coupling between the two nodes. 
It is easily verified that the eigenvalues are 



E±=E± Ve ±l,p = E± y/V 2 - T 2 /4 - iT/2, (6) 

where <f> = arcsin(r/2V). For r — * {<fi ~ » 0) this yields 
the correct eigenvalues E ± V of Hq. Note that for T < 
2V the negative imaginary part of E± is identical for both 



\/2 cos < 



± e ±i<P/2 



(8) 







where the phases <f> depend on T such that in the limit T 
one recovers the eigenstates of Hq, 

We note, however, that finding the bi-orthonormal basis set 
is not necessary for the following calculations. One just has 
to require that the basis sets of H and H^ are orthonormal, 
respectively. In this way one diagonalizes H and H' sep- 
arately, which in the end leads to the same eigenvalues and 
eigenstates as the approach described above. 

When the coupling to the environment vanishes (A — > 0) 
one obtains the survival probability directly from the eigen- 
states and eigenvalues of H . For T < 2V one has 



n(t) 



-ri 



cos 2 (0 + tV cos< 



cos z 



(forA = 0). (9) 



We note that for values T > 2V the dimer is overdamped. 

When considering the dimer without traps (r = 0) but cou- 
pled to the environment, Eq. © simplifies and, from the solu- 
tion for p, one obtains the transition probabilities 
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A sin (WAV 2 - A 2 ) 
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(forT = 0X10) 



and n^l(t) = 1 — n^l (t). For A — > one recovers the sim- 
ple oscillatory behavior of the transition probabilities [namely, 

liniA^o ^i ! (t) = cos 2 (Vt)]. For A > 0, i.e., with coupling to 
the surroundings, the transition probabilities still show oscil- 
lations superimposed on an exponential decay in time which 
tends to the classical equipartition value of 1 /2. 

In order to combine the results of Eq. (O - we will only fo- 
cus on values T < 2V at this stage - and Eq. ([Tol l, we expand 
in both equations all terms except the exponentials to first or- 
der in r and A, respectively. Note that for Eq. (O we obtain a 
product of exp(— Tt) and cos 2 Vt, which is the simple oscil- 
latory behavior of the dimer without trap. Now, the coupling 
to the environment affects all transitions but still conserves 
probabilities. Therefore, we replace the term cos 2 Vt by the 
expansion of Eq. ( fTOl ). by which we obtain 



n(t) 



-r* (o) 
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(11) 
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A. PIMC 



IV. COMPARISON OF LVNE TO PIMC 



In order to corroborate our results we will compare the phe- 
nomenological LvNE approach described above to the numer- 
ically exact PIMC calculations based on a microscopic mod- 
elling of the dissipative environment. This allows to judge for 
which parameter range the approximation by Eq. <j4j delivers 
satisfactory results. 

When coupling the dimer to a bath, the total Hamiltonian 
reads J? to t = H + Hi + Hb, where the dimer-bath cou- 
pling and bath are described in the framework of the Caldeira- 
Leggett model [12], 



£• 



v ( k 



1 



(12) 

Here, P K and X K are the momentum and position operators 
of the bath degrees of freedom, respectively, while m K and 
lo k denote their mass and frequency. The counter term, which 
prevents a renormalization of the free dimer's energy levels 
due to the environmental coupling [111], is absent in Eq. (fT2l 
since it reduces to a physically irrelevant constant in the case 
of a two-level system. After tracing out the environmental de- 
grees of freedom, the onsite population of node \n) becomes 

7r„,i0) = £mS Ht)in exp(~S[a}-$[a} \ , (13) 

where a denotes a closed quantum path in terms of the eigen- 
states of <x z with cr(0) = —1 (referring to the initial prepara- 
tion in node |1)), and S[a] is the action of the free dimer. The 
influence of the environment is summarized in the Feynman- 
Vernon influence functional $[cr] 11311 . which is completely 
determined by the environment's spectral density, 



2 

ifl z — ' rn i,i 



m K uj K 



(14) 



for further details, we refer to Ref . ifTTll . 

As the exact dynamics Eq. ( fL3l can not be calculated ana- 
lytically, one has to resort to a numerical evaluation of the path 
integral. Here, the PIMC method has proven to be a promis- 
ing approach to obtain numerically exact results even in re- 
gions of parameter space where approximative methods fail 
(for details, see e.g. Refs. H d). In our case it is straight- 
forward to adopt the approach presented in Ref. IU5I1 once the 
free dimer's forward and backward propagators, which define 
S[a], are expressed according to 



iE a t/h 



\exp(-iHt/h)\n') = ]T )($>' 

<T = ± 
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er = ± 



(15) 



Figure Q] compares the survival probabilities of a dissipa- 
tive dimer obtained from the approximative LvNE approach 
to the numerically exact PIMC calculations for a bath with 
ohmic spectral density ( TPfl ) with exponential cutoff, J(ui) = 
2-7T 'etui [ e~ ~ u '' !u3a . Here, the onsite energies E and the coupling 
elements V have been taken to be equal, E = V = 1, while 
the temperature is fixed to ksT = V, and we set lo c = 5V. 
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which, even for r ^ 0, are complex conjugate to each other. 



FIG. 1: (Color online) PIMC results (circles) for a dimer with 
r = 0.1 and different system-bath couplings a = X/tt: (a) a — 0, 
(b) a = 1 / 10, and (c) a — 1/4. The solid lines represent the numer- 
ical solution of the LvNE equation, the dashed lines show the cor- 
responding analytical results obtained from Eq. ([9} for a = A = 
and from Eq. jilt for a = 1/10 and a — 1/4. The dotted blue line 
shows the long time limit H(t) ~ exp(— Ft). Panel (d) shows the 
corresponding long-time behavior of the numerical LvNE solution 
for the three different values of a and additionally the behavior for 
large couplings a = 10. 

For small trapping strength (T = 0.1) and vanishing cou- 
pling to the environment (a = 0), Fig.QJa), the PIMC calcula- 
tions coincide with the result of Eq. (|9]). A moderate increase 
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of the coupling (a = 1/10), Fig. Gib), still leads for Eqs. © 
(solid lines) and ( ITTb (dashed lines) to results which are in ex- 
cellent agreement with the findings of the PIMC calculations 
(symbols). When increasing the coupling further to a = 1/4, 
Fig.QJc), however, the approximate solution, Eq. dTTb , begins 
to deviate from the LvNE and the PIMC calculations, which 
are still in very good agreement. 

As the numerical effort of real-time PIMC simulations 
grows exponentially with time, they can cover only short-to- 
intermediate timescales. However, the agreement between the 
LvNE and the PIMC calculations in the weak coupling regime 
allows to compensate for this shortcoming, by using at longer 
times the LvNE results, see Fig.QJd). 

As mentioned earlier, strong couplings A in the LvNE 
lead to the Zeno limit, and therefore disagree with the large- 
coupling/high-temperature behavior of the PIMC formalism. 
Figure Q2d) shows II(i) for a = 10 (dashed-dotted line), 
which clearly deviates from the long-time behavior of the 
curves for small a-values. This corroborates the fact that the 
LvNE in the Lindblad form only yields the correct long-time 
behavior for rather weak couplings to the environment. 




FIG. 2: (Color online) Dimer with a = 1/10 and different trapping 
strengths T. The lines represent the numerical solution of the LvNE 
and the symbols show the corresponding PIMC results. 

For a — 1/10 and varying T, Fig. [2] compares H(t) ob- 
tained from the numerical solution of the LvNE to the results 
of the corresponding PIMC calculations. The trap is part of 
the system, i.e., it directly enters the Hamiltonian H. There- 
fore, unlike a change in a, a change in V leads to no significant 
differences for II (t) obtained from the two methods. Never- 
theless, varying the trapping strength has strong implications 
for the survival probability: At short times (t < 2) the sur- 
vival probability II (t) is smaller for smaller T. At larger times 
the decay is very pronounced for larger T-values; smaller T- 
values do not dampen the oscillations which are superimposed 
on the decay. This effect will be weakened when increasing 
the coupling to the environment. 

V. THE INCOHERENT CASE 

Finally, in the limit of very strong coupling to the environ- 
ment, all coherences will be quickly destroyed and the re- 



sulting transport becomes purely incoherent. The dynamics 
is governed by a master equation whose transfer matrix is 

One obtains the parameters E, V, and T by fitting to the PIMC 
calculations at large temperatures. In principle one can obtain 
the transfer rates V from a golden rule approach 11611 . The 
matrix T has purely real eigenvalues, 

\±=E± Ve^ = E + f/2±\jv 2 + t 2 /A 1 (17) 
Similar to the coherent case, the eigenstates of T read 

i* ±} - v^m (18) 

where ip — arcsinh(r/2V r ). Now, the classical survival prob- 
ability is readily obtained as 

*< fp^r /oi cosh {ib + tV cosh ib) 

P(t) = Pt.i t = e-* (B+r/2) ^ r . ^, (19) 

cosh-0 

which for not too small t gives rise to a simple exponential 
decay with exponent A + . Taking E = V and f « 1 one 
obtains P(t) ~ e~* r / 2 , which is independent of E and V. 
We note the difference in the overall exponential decay of the 
LvNE which was proportional to exp(— Tt). Comparing this 
to the PIMC results for strong couplings a, we obtain V = 2T. 



VI. CONCLUSIONS 

We have studied the trapping of excitations in dimers which 
are coupled to a dissipative environment. By using the (ap- 
proximative) Lindblad form of the LvNE on one hand and the 
(numerically exact) PIMC calculations on the other, we were 
able to specify the range of coupling parameters in which both 
methods agree with each other. Moreover, matching the two 
approaches by appropriately adjusting the coupling parameter 
a in the LvNE allows to extrapolate the short-to-intermediate- 
time PIMC result to - in principle - arbitrary long times. 
Since PIMC is numerically exact, the combination with the 
LvNE is ideal for also studying large systems at long times 
for a broad range of couplings. 
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